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Abstract 

Size-selective concentration of particles in 3D turbulence may be related to collec- 
tion of chondrules and other constituents into primitive bodies in a weakly turbulent 
protoplanetary nebula. In the terrestrial planet region, both the characteristic size 
and narrow size distribution of chondrules are explained, whereas “fluffier” particles 
would be concentrated in lower density, or more intensely turbulent, regions of the 
nebula. The spatial distribution of concentrated particle density obeys multifractal 
scaling, suggesting a close tie to the turbulent cascade process. This scaling behavior 
allows predictions of the concentration probabilities to be made in the protoplanetary 
nebula, which are so large (> 10 3 - 10 4 ) that further studies must be made of the role 
of mass loading. 

Background and Introduction: Most chondritic meteorites are composed in large part of 
mm-sized, once-melted silicate particles (chondrules) and metallic grains which are narrowly 
sorted, most likely by aerodynamic cross section (1); various hypotheses have been advanced 
to explain these properties (2,3). In prior work (3), we noted that chondrule-sized particles 
are unable to settle to the nebula midplane unless turbulence is vanishingly small. Instead, 
we proposed that, following an undefined formation process in the nebula (4), chondrules 
pursue an extended free-floating existence until plausible conditions of nebula gas density 
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and turbulent intensity in the terrestrial planet region concentrate them into much more 
massive, but still not solid, ‘‘clouds" which might have enough coherence to resist disruption 
by eddies and settle to the midplane for subsequent accumulation into planetesimals. We 
think of this process as "primary accretion”. The early stages are easily generalized to a 
wide range of fluffier particles in lower density nebula regions. The final stage (settling of 
dense clouds) remains unstudied, and will require better understanding of the behavior of 
particle ensembles whose density is large enough to affect the gas flow properties. 

We have defined the concentration factor C as the ratio of the local particle density 
to its global average; in numerical studies to date (3, 5, 6), Reynolds numbers of 10 2 — 
10 3 produce C 40 — 300. While C is seen to increase systematically with increasing 
Reynolds number, estimates of concentrations at nebula Reynolds numbers had required 
sizeable extrapolations (3). Here, we utilize our recent discovery of the Reynolds-number- 
independent fractal properties of the particle density field (6) to provide a firmer basis for 
predictions under nebula conditions. 

Turbulence and turbulent concentration (TC) in the nebula: Homogeneous, isotropic, 
3D turbulence is characterized by a cascade of energy through a range of scales, known as 
the inertial range, from the largest (or integral) spatial scale L (with associated velocity 
Vi) to the smallest (or Kolmogorov) scale tj = LRe~ 3 / 4 where it is dissipated (7). The 
intensity of turbulence is often characterized by the Reynolds number Re = LVi,/i/ m , where 
v m is the molecular viscosity; Re is the ratio of transport by macroscopic motions to that 
by molecular motions. This definition expresses the velocity and length scales in terms 
of viscosity; however, TC is based more on turbulent kinetic energy k than on turbulent 
viscosity (8). Therefore, we carefully distinguish between these quantities. We redefine 
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Rt = oc v cH/u m , where c and H are the sound speed and vertical scale height of the nebula, 
and Gr„ is the familiar Shakura-Sunyaev parameter. Similarly, we define cu- by k = V'H 2 = 
\akC z . With the stipulation that true turbulent eddies cannot have frequencies smaller than 
the local orbital frequency, we find L « Hy/ak, and Vl ~ c^/ai t. 

There is current concern about how, or whether, turbulent kinetic energy is maintained 
in neutral nebula gas in the terrestrial planet region (9), which is related to the subtle 
distinction between a„ and a*. If ar„ « a*, for instance, the energy source for k would 
need to be something other than simple viscous evolution of a Keplerian disk. There are 
several energy sources (9) for which typical estimates of a* axe ~ 10 -4 - 10 -2 , and thus 
Re = 10 8 - lO 10 (10). Thus, in spite of pressing concerns, the issue of turbulence appears to 
remain open. Henceforth we neglect the distinction between a* and a„ and treat Re and its 
associated a as providing L and Vt, whether or not their product constitutes a net positive 
Reynolds stress. 

We presume a Kolmogorov- type inertial range, within which each length scale / is charac- 
terized by velocity v(l) = V^l/L) 1 ^ 3 and eddy frequency uj(1) = v(l)/l = Q(l/L)~ 2 ^ 3 , where 
in the nebula ft = <*;(£) .= Vi/ L is the local orbital frequency under solar gravity (3,7). 
Because the most interesting scales for particle concentrations are on the order of q « L, 
deviations from isotropy due to rotation are not a major concern. For particles smaller than 
the geis mean free path (ie, smaller than several cm radius under nebula conditions), the 
stopping time due to gas drag is t , = rp,/cp s , where r and p, are particle radius and internal 
density and p g is the gas density. The particle Stakes number Sti = t,u(l) determines the 
particle response to eddies of a particular scale and frequency; previous studies have shown 
that the optimally concentrated particles have St v = t,u{r]) ss 1 (5), that is, their stopping 
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time is comparable to the Kolmogorov eddy turnover time. 

Our original work (3) concluded that the optimally concentrated particles in the ter- 
restrial planet region of the nebula have radius and density comparable to those of chon- 
drules; furthermore, TC is easily generalized to a wide range of nebula conditions. The 
expression from (3) for optimally concentrated particle size may be rewritten as rp, = 
(mtf 2 /4cr// 2 )°' 5 (I!/a) 0 ' 5 « 1.4 x lO~ 3 R]$j ‘ 5 aT 0S g cm -2 , where I] is surface mass density, 
rn.Hi and are the mass and cross-section of H 2 , and properties of the central scale height 
in a Hayashi-type “minimum mass” nebula are assumed (11). This relationship is shown in 
figure 1 for several different typical locations. Note that porous aggregates (PA), having 
considerably lower radius-density product than chondrules, are optimally concentrated at 
the low gas densities which characterize the outer planet region or low density regions high 
above the nebula midplane. Such porous objects are easily produced because the low relative 
velocities in turbulence of low-density grains and their constituent monomers lead to large 
sticking efficiency and minimal disruption (13). TC may thus have been ubiquitous, initi- 
ating the formation of “cometesimals” from porous grain aggregates at 10-30 AU. However, 
textural evidence might be difficult to obtain from this regime; subsequent compaction would 
obliterate evidence for any preferred size or density of easily squashed fluffy constituents. 
Chondrules and their parent meteorites, by nature of their availability and unique, persistent 
textures, provide the most obvious initial testing ground for turbulent concentration. 

Concentrated particle size distribution: 

We have recently determined the properties of turbulent concentration from direct numer- 
ical simulations of homogeneous, isotropic, incompressible 3D turbulence at three Reynolds 
numbers Re = 110, 426, and 1300 (6,14). The particles can be given arbitrary aerodynamic 
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Preferred radius-density product r*rho vs 'alpha' 



nebula "alpha" 

Figure 1: The optimal Stokes number for turbulent concentration, converted to the product 
of particle radius and density by adopting nominal nebula gas density and local rotation 
rate at various locations in a typical solar nebula, as functions of a. Solid silicate particles 
with chondrule sizes (CH) are favored in the terrestrial planet region for nominal a (3). 
Smaller gas densities at larger distances from the sun concentrate smaller rp p products 
(porous aggregates, or PA) for any a. It has also been suggested that extremely low density, 
intensely turbulent regions may concentrate tiny grains or very fluffy aggregates (FA) (12). 
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stopping times their motions respond only to gas drag and are integrated in the spa- 
tial domain. The computationally intensive calculations are run on 16 C'90 cpus at Ames 
Research Center. 

One important finding relates to the detailed form of the size distribution of selectively 
concentrated particles. If dense clouds of particles are precursors of primitive bodies, the 
distribution with Stokes number of particles in dense regions should correspond to the dis- 
tributions found in chondrites. Simulations were initiated with uniform spatial distributions 
of particles, themselves uniformly distributed in stopping time over the range St v = 0.1-6 
(14), and the relative equilibrium abundance of particles was studied as a function of Stj, 
and C. In the large-C limit of interest, the shape of the distribution of relative abundance vs 
St v was found to be independent of both C and Re (14). Thus, the numerical results should 
be valid as a prediction of the size distribution in clouds under nebula conditions. 

In figure 2 we compare a histogram of relative particle numbers as a function of St 
from our numerical simulations, with binned data for 253 disaggregated chondrules from 
the L4 ordinary chondrite (OC) ALH85033. While a few ostensibly similar data sets exist 
in the literature (15), data for this and several other meteorites imply that merely measur- 
ing the chondrule radius distribution and assuming some mean density misrepresents the 
rp a distribution because of chondrule- to-chondrule density variations (16). The meteorite 
histogram was aligned horizontally with the predicted histogram by assuming that the opti- 
mally concentrated particle size lies at St v = 1. The shape of the theoretical profile contains 
no adjustable parameters, and the data and model profiles are in excellent agreement. The 
implication is that TC by itself can explain the very narrow chondrule size distribution, what- 
ever the chondrule formation process may have produced. Larger or smaller “chondrules” 
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Figure 2: Comparison of the size-density distribution of 253 chondrules disaggregated from 
L4 ALH 85033 (solid symbols) with the theoretically determined shape distribution (open 
symbols) for the relative abundance of particles concentrated in turbulence (14,16). The 
theoretical profile is independent of C and Re, close to lognormal in shape, and a good 
match to the shape of the data. Four different meteorites now exhibit agreement this good 
(Cuzzi et al 1999) 
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might have existed, but would simply not be concentrated to such a significant extent for 
delivery to the midplane unless swept up into chondrule rims or broken to the proper size. 
Of course, some particles of arbitrary size will always appear just by accident, and parent 
body processes will confuse the situation further. 

Particle concentration is a multifractal with Re -independent properties: A nonspecjalist’s 
definition of a fractal is a structure which is generated by sequential application of a rule on 
regularly decreasing spatial scales. Simple fractals with constant (but non-integer) dimension 
result from rules which produce a binary distribution of the local density (say, either 0 or 1) - 
and are invariant to changes in scale. Examples of these are the Cantor set or the Sierpinski 
gasket, in which segments of a line, or portions of a plane (17), are simply removed without 
changing the surrounding values. Their average density, as a function of scale e, may be 
written as p(e) — poi~ D with dimension D = 0.63 and D = 1.52 respectively. 

For comparison, multifractals result from application of rules in which the local measure 
is changed while conserving the total measure - for example, by unequal (but invariant) 
repartitioning of the content of a bin into sub-bins at each of the many subsequent stages 
of a cascade process (17). Such quantities have no well defined local value in the limit 
of diminishing bin size; that is, their local values are spatially spiky or “singular” (17). 
Considerable study has been devoted to multifractals in turbulence, because turbulence 
(and specifically its inertial range) is the archetype of a cascade process (7,17). For instance, 
dissipation of turbulent kinetic energy, which occurs on the Kolmogorov scale, is not spatially 
uniform but has the spatial distribution of a mutifractal (18). The spatial distributions of 
multifractals are predictable in a statistical sense, using probability distribution functions 
(PDFs) which are derived directly from their dimensions. 
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In multifractals, the dimension varies with the value of the measure: in our case, the 
measure is particle concentration C which is defined in terms of e“, where a is a scaling 
index which can be regarded a s a local dimension for C. The key statistical descriptor of 
a multifractal is its singularity spectrum f(a), where f(a) plays the role of a dimension of 
the probability distribution function for C (19,18,6) and is the key element in the PDF’s 
(discussed further below). 

To make predictions at values of Re very different from directly accessible values, thus, 
the Re— dependence of f(a) is crucial. Using the methodology of (18), we have shown that 
the spatial distribution of optimally concentrated (St, = 1) particles is a multifractal, with 
f(a) invariant over more than an order of magnitude in Re (i?e=100 - 1300) (6). This 
/?e-invariance for f(a) is only seen if binning is tied to some fundamental flow-relative scale 
such as the Kolmogorov scale tj, which varies with Re in a known way. Dissipation has also 
been shown to have an Re- independent /(a) - from numerical work at Re ~ 100, including 
our own, through laboratory experiments with Re ~ 10 4 , to experimental studies of the 
atmospheric boundary layer with Re ~ 10 6 (20,18). Particle concentration and dissipation 
are physically connected through their preference for the Kolmogorov scale, and the shapes 
of the singularity spectra for dissipation, and for particle concentration binned at close to 
the Kolmogorov scale, are similar (6). These facts suggest that the physics which leads to 
turbulent concentration is closely related to the turbulent cascade process - as for dissipation. 
The turbulent cascade is known to have statistically predictable, /^-independent properties 
in the inertial range. Based on these arguments, we believe and presume that the particle 
concentration singularity spectrum f(a) remains /2e-independent to far larger values than 
those of our numerical experiments. Given this invariance, we can predict nebula conditions 
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more confidently than from extrapolation alone, as was done in (3). 

The scaling index a is related to C as follows: the fractional probability P, of a particle 
lying in a bin which contains ;V, particles out of N p total particles is defined as P, = NJ N p = 
e a ‘, where t = bin size/domain size. The associated concentration factor C, = - where u, 

is the volume of a bin and v is the total volume of the domain; thus C, = — (e“'/ e3 ) = 

e a, ~ 3 . As mentioned above, Re — invariance of f(a) is only expected when the bin size is some 
multiple J of ij. The domain normalized bin size is then e = Jr//ML = -fa Re ~ 3/M = 1 /7£, 
where ML is the domain size (18) expressed in units of the integral scale Z-, and we use the 
inertial range relationship tj = LRe~ 3 ^ 4 . Thus, C = 7 Z 3 ~ a , or a = 3 — For a nebula 
characterized by some turbulent a, with bin size 2tj, integral scale L = H y/a, and a domain 
size of 2 H (21), one obtains 

K = a l/4 {cHlv m f ,A . (1) 

The normalized PDF for a is usually written a s F(a ) = p(a)e 3 ~P a \ where the prefactor 
p{a) can be approximated as \J\n7 Z (18); We define the PDF F V {C) as the volume fraction 
occupied by bins having concentration factor C, with /c"“ F v (C)dC = f*™” F (a) da = 1. 
Transforming variables and their PDFs, we get 


F,(C) = F(a) | ^ |= 




( 2 ) 


dC ' Clnft''' CV WR 

Note how the function f(a) assumes the role of a dimension for F V (C). The fraction 
of particles occupying bins with C is F P (C) = CF(C), and both F V (C ) and F P (C ) have 
cumulative versions F v (> C) — F v (C)dC, and F p (> C) = F p (C)dC. We have ob- 
tained numerical results which validate the (ergodic) assumption that the fraction F p (> C ) 
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(spatially averaged over all particles at several snapshots in time) is equal to the fraction of 
time F t {> C) spent by a given particle in regions denser than C (temporally averaged over 
extended trajectories for a few particles). 

If /(a) is indeed a Re- independent, universal function for optimally concentrated (St n = 
1) particles, we can then predict PDF’s for any given nebula a(= a*) using equations (1) 
and (2) derived above. As a check on the method (there are subtle normalization issues), an 
average /(a), obtained from our numerical calculations at ail three Re values (22), was used 
to calculate F p (> C ) (= F t (> C)) for comparison with the distributions directly determined 
from numerical results at our three Re values (6). As shown in figure 3, the single /(a) does 
quite well at predicting the PDFs at all three Re. Also shown in figure 3 are associated 
predictions of F p (> C ) for four values of nebula a, using equations (1) and (2). Based on 
these predictions, St = 1 particles spend 10 percent of their time in regions with C > 10 3 , 
and 1 percent of their time in regions with C > 10 4 , under nebula conditions. These results 
might help us understand why some chondrule types seem to have been formed under highly 
oxidizing conditions (1). 

Encounters with dense clouds, and entrapment: The PDF’s can then be combined with 
the particle velocity through space V p to calculate the “encounter time” T enc of a chondrule 
with a region of arbitrary C, using a duty cycle argument: F t (> C ) , 

where £,„(> C) « T enc is the time spent traversing a bin with concentration greater than 
C. Thus, for bins of dimension 2t) and particle velocity V p , tin (> C) = 2rj/V p , SO Tenc « 

= VrF^Fc)' 411(1 the encounter rate is it: = ^ F p(> C )• We have verified 
numerically that, as expected for particles with stopping times f, much shorter than the 
overturn time of the largest eddies fi _1 (23,11), V p is nearly identical to the typical gas 
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Concentration factor C 

Figure 3: Probability Distribution Functions (PDF’s) for the fraction of particles lying in 
regions with concentration factor greater than C ( F p (> C)), or, equivalently, the fraction 
of time spent by any particle in such regions (i*t(> C )) as discussed in the text. The three 
sets of points are binned directly from our numerical simulations; the associated curves are 
calculated from a single averaged f(a) obtained from all three values of Re. The curves 
without points use the same f(a) to predict PDF’s at the larger Re corresponding to four 
plausible nebula a values: 10 -2 (solid line), 10 -3 (long dashed line), 10 -4 (short dashed line), 
and 10 -5 (dotted line). 
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velocity V 3 . where V 0 w ^ \f2k % c^/cv (23). 

Specifically, we may calculate the encounter rate (and time) with a cloud so dense that 
a particle becomes entrapped with its neighbors and removed from further free circulation. 
Normally, particles traverse dense regions without incident, as their gas drag stopping time 
is longer than their transit time (3). An entrapment threshold occurs when interparticle 
collisions prevent particles from passing through a cloud; this implies a critical cloud optical 
depth Tcou of unity, defining a critical C co it- For small, dense 2»7-sized clumps, t co u = 1 = 
‘2vCcou(Pch/m)*r 2 , where p c h is the average (unconcentrated) chondrule mass density, and 
m and r are chondrule mass and radius. Then C co n ~ 300 rp 3 /Tjp g f ch m 3 x 10 5 / /«.*, where 
f c h is the fraction of all solid mass in chondrule form, p 3 ~ 3 is the density of a chondrule, 
and nominal parameters at 2.2 AU are assumed (11). Such clouds are orders of magnitude 
denser than the gas ( Pd/Pg ~ 1500 for / c * = 0.3 and C co ii ~ 10 6 ). For these parameters, 
and the PDF’s we have in hand, T enc is on the order of 10 4 years - short compared to nebula 
evolution timescales. This quantifies the first part of our earlier suggestion (3) that such 
an encounter at T enc (Ccoii ) ends the freely wandering life of a chondrule. A second stage, 
in which the chondrule and its cluster of neighbors might settle toward the midplane under 
the vertical component of solar gravity as a dense, but still far from solid, unit (3), remains 
more speculative and is under study. Alternately, TC might merely provide a large increase 
in the collisional aggregation rate of optimally sized particles. 

Discussion: 

Limitations due to mass loading: Even with this approach being clear in principle, there 
is one important unknown factor remaining. The cascade process model is only valid as 
long as no “new” physics emerges at some step in the cascade to change the “rule”, or 
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f[n), applicable to subsequent steps in the cascade. However, the above predictions imply 
that particle collisions trap particles in small, dense clouds only when p p / p g > 10 3 . While 
these high concentrations literally relate to regions comparable to tj in size, where there 
is no turbulence to damp, the cascade that produced such a cloud from its "penumbra' 1 , 
must have extended to larger sizes within which turbulent motions may be damped by 
lower concentrations. We have made some preliminary calculations of f(a) for mass loading 
Ppl Pg 1. Turbulent concentration persists, but its f(a) is altered in the sense that high C 

values have a lower probability. Clearly, this effect must be quantified before more specific 
predictions of T enc and other accumulation timescales can be made. 

Complications due to cloud evolution and parent body processing: The subsequent evolu- 
tion of dense clouds in the vertical component of solar gravity must still be studied. Do they 
retain their identity or disperse? In subsequent stages, solids must be compacted to orders 
of magnitude higher density than the dense clouds produced by turbulent concentration - 
possibly due initially to collisions between dense clouds en route to, or within the midplane, 
followed by midplane sweepup of the dense clouds or their constituents by objects which can 
begin as small as ten meters (ie (11), section 5.4; also (13)). The fact that most primitive 
meteorites contain ample evidence for abrasion, fragmentation, and other mechanical pro- 
cesses which may well have continued long after the nebula gas vanished and aerodynamic 
processes became irrelevant, warns us that we should not expect turbulent concentration to 
explain ail properties of even “primitive” meteorites. 
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because we are interested in low-probability events; the distinction is discussed by Chhabra, 
A. B. and K. R. Sreenivasan (1991) Phys. Rev. A. 43, 1114-1117. Our /(a) is well fit by 
the function /(a) = —12.414 + 89.659/a — 132.01/a 2 . 

23) Volk, H. J., F. C. Jones, G. E. Morfill, and S. Roser (1980) Astron. Astrophys. 85, 
316; also Markievicz, W. J., H. Mizuno, and H. J. Volk (1991) Astron. Astrophys. 242, 
286-289 
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